#!/usr/bin/octave

sizes = [35, 50, 75, 100, 139, 150, 160, 200, 250, 300, 340, 350, 400, 450, 490, 500];

arg = argv;
standard = load("-ascii", arg{1})';
sample = load("-ascii", arg{2})';

% pick the highest peak in the last two-thirds of the run
height = max(standard(floor(end/3):end)) - min(standard)
quantum = height/10
threshold = 2*quantum
x = 1:length(standard);
peak = find(standard > threshold);
first = floor(sum(peak(1:2))/2) - 100 # The first peak is garbage; this used to be peak(1) - 200
lastpeak = peak(end)
last = lastpeak + 150; % don't want to mangle the last peak

% clip the high levels in the leading garbage, and clip the low levels everywhere
outliers = find(standard > 1.5*height);
standard(outliers) = height;
outliers = find(standard < 0.5*threshold);
standard(outliers) = 0;

lastpeak = lastpeak - first;
scale = max(sizes)/lastpeak;
standard = standard(first:last);
sample = sample(first:last);
x = (1:length(standard)) * scale;

ref = zeros(1, length(standard));
for i = (1 : length(sizes))
  ref = ref + height * pulse ( sizes(i), 0.18, x ); 
endfor

#from = find(x > 0)(1);
#to = find(x > 120)(1);
from = 1;
to = length(x);
#sample = 1:length(standard);
X = vertcat(x, ref, standard, sample);
X = X(:, from:to)';
save("-ascii", sprintf("%s", arg{3}), "X");

